
#ifndef ROPUFU_MULTI_SENSOR_STATISTIC_PARALLEL_MAX_CUSUM_HPP_INCLUDED
#define ROPUFU_MULTI_SENSOR_STATISTIC_PARALLEL_MAX_CUSUM_HPP_INCLUDED

#include "../information/information_none.hpp"
#include "../more_math.hpp"
#include "../window/window_fixed.hpp"
#include "../window/window_none.hpp"
#include "cusum.hpp"

#include <array>
#include <cmath>
#include <cstdint>
#include <ostream>
#include <string>
#include <type_traits>
#include <vector>

namespace ropufu
{
    namespace statistic
    {
        namespace detail
        {
            template <typename t_window_type>
            struct parallel_max_cusum_default_window
            {
                static t_window_type make()
                {
                    return t_window_type();
                }
            };

            template <std::size_t t_count_streams>
            struct parallel_max_cusum_default_window<window_fixed<t_count_streams>>
            {
                static window_fixed<t_count_streams> make()
                {
                    return window_fixed<t_count_streams>(1);
                }
            };
        }

        // WARNING: Only implemented for independent streams.
        template <std::size_t t_count_streams, typename t_information_type = void, typename t_window_type = window_fixed<t_count_streams>>
        struct parallel_max_cusum
        {
            static const std::size_t count_streams = t_count_streams;
            typedef typename switch_void<t_information_type, information_none<count_streams>>::type information_type;
            typedef typename switch_void<t_window_type, window_none<count_streams>>::type window_type;
            typedef parallel_max_cusum<count_streams, t_information_type, t_window_type> type;

        private:
            information_type m_prior_knowledge;
            window_type m_window;
            std::vector<cusum> m_local_cusum;
            double m_current_statistic = 0.0;
            bool m_has_reset = false;

            parallel_max_cusum() = delete;

        public:
            explicit parallel_max_cusum(const information_type& prior_knowledge)
                : m_prior_knowledge(prior_knowledge), m_window(detail::parallel_max_cusum_default_window<window_type>::make()), m_local_cusum(prior_knowledge.cardinality())
            {
            }

            explicit parallel_max_cusum(const information_type& prior_knowledge, const window_type& window)
                : m_prior_knowledge(prior_knowledge), m_window(window), m_local_cusum(prior_knowledge.cardinality())
            {
            }

            bool has_reset() const
            {
                return this->m_has_reset;
            }

            void reset()
            {
                this->m_window.reset();
                for (cusum& local_stat : this->m_local_cusum) local_stat.reset();
                this->m_current_statistic = 0.0;
                this->m_has_reset = false;
            }

            double observe(const std::array<double, count_streams>& instantaneous_loglikelihood)
            {
                this->m_current_statistic = 0.0;
                this->m_window.observe(instantaneous_loglikelihood);
                std::size_t sz = this->m_window.get_depth();

                this->m_has_reset = sz == 0; // TODO: implement as <window_type> property.
                if (sz == 0) // If it's a reset point...
                {
                    for (cusum& local_stat : this->m_local_cusum) local_stat.reset();
                    return 0.0;
                }

                // TODO: implement non-unit weights (right now p_{\cA} = p^{|\cA |} / \sum p^{|\cB |} with p = 1).
                double log_weight_norm = std::log(static_cast<double>(this->m_prior_knowledge.cardinality()));
                std::size_t index = 0;
                for (const std::array<bool, count_streams>& subset_mask : this->m_prior_knowledge)
                {
                    double subset_loglikelihood = 0.0;
                    for (std::size_t k = 0; k < count_streams; k++)
                    {
                        // WARNING: This is where independence kicks in.
                        if (subset_mask[k]) subset_loglikelihood += instantaneous_loglikelihood[k];
                    }
                    double local_result = this->m_local_cusum[index].observe(subset_loglikelihood) - log_weight_norm;

                    if (index == 0) this->m_current_statistic = local_result;
                    else if (local_result > this->m_current_statistic) this->m_current_statistic = local_result;
                    index++;
                }
                return this->m_current_statistic;
            }

            // Returns the complete loglikelihood ratio.
            const std::array<double, count_streams>& complete() const
            {
                return this->m_window.complete();
            }

            // Returns the accumulated loglikelihood ratio.
            const std::array<double, count_streams>& accumulated() const
            {
                return this->m_window.accumulated();
            }

            friend std::ostream& operator <<(std::ostream& os, const type& that)
            {
                return os << "parallel_max_cusum_" << that.m_prior_knowledge << "_" << that.m_window;
            }
        };
    }
}

#endif // ROPUFU_MULTI_SENSOR_STATISTIC_PARALLEL_MAX_CUSUM_HPP_INCLUDED
